## Bar plots decomposing fiscal income and expenditure in the GCC either since 2010 or for certain snapshots

### Housekeeping: setting up R environment
knitr::opts_chunk$set(echo = TRUE, warning = FALSE)
Sys.setenv(LANG = "en")
rm(list=ls())
requiredPackages <- c("tidyverse", "ggplot2", # base packages
                      "haven", # read Stata files
                      "readxl", # reading Excel files
                      "viridis", # colors
                      "scales", # for color palettes
                      "ggpattern" # for pattern-type plots
) 

ipak <- function(pkg){
  new.pkg <- pkg[!(pkg %in% installed.packages()[, "Package"])]
  if (length(new.pkg))
    install.packages(new.pkg, dependencies = TRUE)
  sapply(pkg, require, character.only = TRUE)
}
options(scipen=999)
ipak(requiredPackages)
rm(list=ls())

## Barplot decomposing yearly fiscal revenue and expenditure in GCC
Bahrain <- read_dta("../DATADIR/Bahrain_raw.dta")
Bahrain$CountryName <- "Bahrain"
Bahrain$ISO <- "BHR"
Kuwait <- read_dta("../DATADIR/Kuwait_raw.dta")
Kuwait$CountryName <- "Kuwait"
Kuwait$ISO <- "KUW"
Oman <- read_dta("../DATADIR/Oman_raw.dta")
Oman$CountryName <- "Oman"
Oman$ISO <- "OMN"
Qatar <- read_dta("../DATADIR/Qatar_raw.dta")
Qatar$CountryName <- "Qatar"
Qatar$ISO <- "QAT"
SaudiArabia <- read_dta("../DATADIR/SaudiArabia_raw.dta")
SaudiArabia$CountryName <- "SaudiArabia"
SaudiArabia$ISO <- "KSA"
UAE <- read_dta("../DATADIR/UAE_raw.dta")
UAE$CountryName <- "UAE"
UAE$ISO <- "AED"

full_dta_q <- bind_rows(Bahrain, Kuwait, Oman, Qatar, SaudiArabia, UAE)

full_dta_q$year <- format(as.Date(full_dta_q$Date, format="%Y-%m-%d"),"%Y")

# Normalize the data by deflating and expressing in USD
exchange_rates_2018 <- full_dta_q %>%
  filter(Date == as.Date("2018-12-01")) %>%  
  select(CountryName, USD_forex_rate_2018 = USD_forex_rate)  

# Join to the main dataset
full_dta <- full_dta_q %>%
  left_join(exchange_rates_2018, by = "CountryName") %>%  
  group_by(year, ISO, CountryName) %>%
  summarise(
    y_ttl_nominal = sum(y),
    g = sum(g / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6)),
    t = sum(t / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6)),
    tax_vat = sum(tax_vat / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6)),
    tax_oil = sum(tax_oil / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6)),
    inv_govt = sum(inv_govt / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6)),
    g_spend = sum(g_spend / GDPdeflatorNSA * 100 / USD_forex_rate_2018 / (population/1e6))
  ) 


full_dta <- full_dta %>%
  dplyr::mutate(
    g_spend = ifelse(is.na(inv_govt), NA, g_spend),
    g = ifelse(!is.na(inv_govt), NA, g),
    t = ifelse(!is.na(tax_oil), NA, t)
  )
  


df_2023 <- full_dta %>%
  filter(year == 2023)

df_2023$year <- NULL
df_plot <- df_2023 %>%
  pivot_longer(
    cols = c("g","t","tax_vat","tax_oil","g_spend","inv_govt"),           
    names_to = "var",              
    values_to = "value"           
  ) %>%
  drop_na() %>%
  # Create a variable indicating whether govt. budget can be split
  dplyr::mutate(fallback = ifelse(var == "g" | var =="t", "stripe", "none"),
                g_or_t = ifelse(var == "g" | var == "g_spend" | var == "inv_govt", "g", "t"))


color_map <- c(
  "tax_vat"           = "#6bded6", 
  "t"           = "#4692C9",  
  "tax_oil"  = "#2171b5", 
  
  "inv_govt"      = "#febd3c",  
  "g"       = "#F0833A",  
  "g_spend"= "#e6550d"   
)


# Re-order 'var' for plotting
df_plot$var <- factor(df_plot$var, levels = c("g", "t", "inv_govt", "g_spend", "tax_vat", "tax_oil"))

# Replace AED with UAE
df_plot$ISO <- ifelse(df_plot$ISO == "AED", "UAE", df_plot$ISO)

final_plot <- ggplot(df_plot) +
  geom_bar( aes(x = g_or_t, y = value, fill = var), stat = "identity") +
  scale_pattern_manual("Position", values = c("none","stripe", "none")) + 
  facet_wrap(~CountryName, nrow=1) +
  #scale_y_log10() +  # Apply log scale to y-axis
  scale_fill_manual(values = color_map,
                    name = "Tax & Spending Components",
                    labels = c("g" = "Total Govt. Exp.",
                               "t" = "Total Taxes", 
                               "inv_govt" = "Govt. Investment", 
                               "g_spend" = "Govt. Spending", 
                               "tax_vat" = "Tax Income", 
                               "tax_oil" = "Oil Revenue"
                                )) +
  labs(
    x = "Country",
    y = "Value (in real 2017 USD)",
    fill = "Segment",
    title = "Sovereign Income & Expenditure per capita in 2023"
  ) +
  theme_minimal() +
  theme(
    legend.position = "bottom",  
    legend.text = element_text(size = 9),  
    legend.key.size = unit(1, "cm"), 
    legend.spacing.y = unit(0.5, "cm"), 
    axis.text.x = element_blank(),  
    strip.text = element_text(face = "bold"),  
    plot.title = element_text(size=22, hjust = 0.5)
  )

final_plot

# Export
ggsave("../RESULTS/govbudget_decomp.pdf", final_plot, 
       width = 8, height = 5, dpi = 300, device = cairo_pdf)




#### A bar decomposition plot of the annual sovereign for Oman
Oman_StatData_PublicFinance <- read_excel("../ORIGDIR/Novel_data/Oman/Oman_StatData_PublicFinance.xlsx", 
                                          sheet = "Tabelle1")

# Define new grouped income and expenditure components
income_total <- "Total_Revenue"
income_groups <- c("Non_Oil_Revenues")  

expenditure_total <- "Total_Expenditure"
expenditure_groups <- c("Defence_&_Security_Expenditures", "Interest_on_Loans", 
                        "Investment_Expenditure", "Civil_Ministries_(Total)",
                        "Subsidies_to_Public_Authority")  

# Standardize column names: Replace spaces with underscores
Oman_StatData_PublicFinance <- Oman_StatData_PublicFinance %>%
  rename_with(~ str_replace_all(., " ", "_"))  

# Adjust for GDP deflator 
Oman_StatData_Quarterly <- Oman_StatData_PublicFinance %>%
  mutate(across(
    c(income_total, expenditure_total, income_groups, expenditure_groups, 
      Net_Oil_Revenues, Natural_Gas_Revenues),  
    ~ .x * GDPdeflator, 
    .names = "adj_{.col}"  
  ))

# Aggregate from quarterly to annual
Oman_StatData_Annual <- Oman_StatData_Quarterly %>%
  mutate(Year = lubridate::year(Date)) %>%
  group_by(Year) %>%
  summarise(
    across(starts_with("adj_"), sum, na.rm = TRUE),  
    
    # Merge Oil Revenues
    Oil_Revenues = sum(adj_Net_Oil_Revenues, na.rm = TRUE),
    
    # Merge Gas Revenues
    Gas_Revenues = sum(adj_Natural_Gas_Revenues, na.rm = TRUE)
  ) %>%
  rename_with(~ gsub("adj_", "", .), starts_with("adj_"))  

# Function to prepare grouped data
prepare_decomp_data <- function(df, total_var, part_vars, group_label) {
  part_data <- df %>%
    select(Year, all_of(total_var), all_of(part_vars)) %>%
    pivot_longer(cols = all_of(part_vars), names_to = "Component", values_to = "Value") %>%
    mutate(Group = group_label)
  
  discrepancy_data <- df %>%
    select(Year, all_of(total_var)) %>%
    mutate(
      Component = "Other",
      Value = df[[total_var]] - rowSums(df[part_vars], na.rm = TRUE),
      Group = group_label
    ) %>%
    select(Year, Component, Value, Group)
  
  # Combine
  bind_rows(part_data, discrepancy_data) %>%
    select(Year, Component, Value, Group)
}

# Prepare Income and Expenditure Data
income_data <- prepare_decomp_data(Oman_StatData_Annual, income_total, c(income_groups, "Oil_Revenues", "Gas_Revenues"), "Income")
expenditure_data <- prepare_decomp_data(Oman_StatData_Annual, expenditure_total, expenditure_groups, "Expenditure")

# Combine datasets
budget_data <- bind_rows(income_data, expenditure_data)

# Ensure Component is a factor to maintain color mapping
budget_data$Component <- factor(budget_data$Component, levels = unique(budget_data$Component))

# Generate a distinct color palette
num_components <- length(unique(budget_data$Component))
custom_colors <- viridis::viridis(num_components, option = "plasma")


budget_data$Component_label <- str_replace_all(budget_data$Component, "_", " ")

# Plot the budget decomposition
Oman_budget <- ggplot(budget_data, aes(x = Group, y = Value/1000000, fill = Component_label)) +  
  geom_bar(stat = "identity", position = "stack") +  
  facet_wrap(~Year, nrow = 1) + 
  scale_fill_manual(values = custom_colors, name = "Component") + 
  labs(
    title = "Oman: Sovereign Budget Decomposition (Annual, 2014-2019)",
    x = NULL,
    y = "Real (2018) Million OMR"
  ) +
  theme_minimal() +
  theme(
    strip.text = element_text(face = "bold", size = 12),
    axis.text.x = element_text(size = 10),
    plot.title = element_text(size = 20, hjust = 0.5),
    legend.position = "bottom"
  )

# Show plot
Oman_budget

# Export
ggsave("../RESULTS/Oman_budget.pdf", Oman_budget, 
       width = 11, height = 7, dpi = 300, device = cairo_pdf)




